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Abstract 



ri ' We describe a further development of the stochastic state selection 

method, which is a kind of Monte Carlo method we have proposed in 
"j^ . order to numerically study large quantum spin systems. In the stochastic 

"*r? ' state selection method we make a sampling which is simultaneous for 

many states. This feature enables us to modify the method so that a 
j^ ■ number of given constraints are satisfied in each sampling. In this paper 

we discuss this modified stochastic state selection method that will be 

called the constrained stochastic state selection method in distinction from 

'■^ ' the previously proposed one (the conventional stochastic state selection 

^ ' method) in this paper. We argue that in virtue of the constrained sampling 

O . some quantities obtained in each sampling become more reliable, i. e. their 

I ^1 ' statistical fluctuations are less than those from the conventional stochastic 

state selection method. 
^SJ . In numerical calculations of the spin-1/2 quantum Heisenberg anti- 

^ ' ferromagnet on a 36-site triangular lattice we explicitly show that data 

0^ , errors in our estimation of the ground state energy are reduced. Then we 

r*^ ■ successfully evaluate several low-lying energy eigenvalues of the model on 

^^ [ a 48-site lattice. Our results support that this system can be described 

by the theory based on the spontaneous symmetry breaking in the semi- 
t~^ ' classical Neel ordered antiferromagnet. 

o 

oo 

O ■ 1 Introduction 

- It is widely recognized that numerical methods based on the first principle are 

r> I quite important in the study of quantum spin systems. Actually the quantum 

S ■ Monte Carlo method has contributed toward enlarging our knowledge of non- 

frustrated quantum spin systems, especially of the spin-1/2 quantum Heisenberg 
antiferromagnet on bipartite lattices [U [51 [3] . Yet it is well known that this Monte 
Carlo method has a limited power in calculations of two-dimensional frustrated 
spin systems. In this situation many studies on numerical methods have been 
proposed and developed, sometimes depending on approximations. One of these 
studies is the coupled-cluster method [H [5], which is a kind of variational method. 
Another approach is the stochastic reconfiguration method[5] whose origin is 



assigned to the fixed- node Monte Carlo method[71[S]. One should also attend to 
various studies on the density matrix renormalization group (DMRG) method[9l 
[TUl [TT] which extend the original DMRG method on a chain [T^l US] to higher 
dimensional systems. The path-integral renormalization group method[T3] is 
interesting as well. 

Recently we have developed another Monte Carlo method, which we call the 
stochastic state selection (SSS) method[15l[16l[17l[18l[T9l[20]. This method has 
a good property common to the ordinary Monte Carlo method that in principle 
one does not need any approximations specific to the system we investigate. 
The sampling algorithm is, however, quite different from the ordinary one since 
the SSS method is based on not importance sampling but a new type of stochas- 
tic selection. In the SSS method we use an operator which generates sampled 
states from any given state. This operator includes a set of stochastic variables 
which are as many as the number of basis states of the vector space under con- 
sideration. The essential point of the selection is that many of these stochastic 
variables are valued to be zero while their statistical averages are all equal to 
one. Therefore in this algorithm we select a relatively small number of basis 
states from a vast vector space in a mathematically justified manner so that the 
statistical averaging processes give us the correct value of any inner product. 
So far all these stochastic variables have been generated independently. 

In this paper we propose the constrained SSS method — a modified SSS 
method whose samples satisfy a number of given constraints. In order to make 
it possible we introduce some dependencies among the stochastic variables stated 
above. From theoretical point of view all we need to restore the original state 
from sampled states is that the statistical average of each stochastic variable is 
equal to one. Therefore we can define some variables in the set as functions of 
other variables, instead of requesting they should be independently generated, 
if each of these functions assures that the statistical average of the dependent 
variable equals to one. Keeping this in mind we develop the constrained SSS 
method, where values of special inner products to represent constraints are 
unchanged in each sampling. We then argue that sampling errors in numerical 
studies can be reduced with suitable constraints and with suitable dependent 
variables. 

As a concrete example, we calculate energy eigenvalues around the ground 
state of the spin-1/2 quantum Heisenberg antiferromagnet on an TV^-site trian- 
gular lattice. The Hamiltonian of the system is 

H^JaJ2s,-Sj, (1) 

where Si denotes the spin-1/2 operator on the i-th site and the sum runs over 
all Ni,{— 3Ns) bonds of the lattice. The coupling J^ is set to 1 throughout this 
paper. In our calculations we employ the power method. The constrained SSS 
method is used to calculate the expectation values of powers of the operator Q, 

Q = li~H, (2) 



where / denotes the identity operator and Z is a positive number which depends 
on the lattice size. For each lattice size we choose one value of I which ensures 
that the ground state eigenvalue of the system corresponds to the eigenvalue of 
Q whose absolute value is the largest. Detailed explanations for it are given in 
Appendix A. 

One reason why we study this system here is that, as is well-known, it is 
a typical example of strongly frustrated systems in two dimensions. Another 
reason is that there has been a long history of investigations into what state 
is realized on the triangular lattice. Lots of studies on this system[31 HI \5\ [TTJ 
imiSlISlIllIinilinillHSHlllSllSni indicate the ground state with the three- 
sublattice order pi]. One should, however, also note recent works [511 [5^ [55] 
which suggest that this quantum system has a richer phase structure than the 
one expected from the classical spin wave theory, as well as other studies [4l [30] 
which show that this system is near the quantum critical point. 

The plan of this paper is as follows. In the next section we describe the 
method. Subsection 12.11 is devoted to a brief review of the conventional SSS 
method [13 . In subsections 12.21 and 12.31 we make a detailed account of the con- 
strained SSS method. First we give a simple example to explain how we impose 
a constraint in a sampling in subsection 12.21 Then extensions to more general 
cases are discussed in subsection 12. 31 Sections [51 H] and [5] are for applications of 
the method to the spin 1/2 quantum Heisenberg antiferromagnet on triangular 
lattices. In section [5] we study the model on a 16-site lattice, where the exact 
eigenstate can be easily obtained. Using this exact eigenstate we evaluate expec- 
tation values of the m-th power of the operator Q by means of the constrained 
SSS method. We find that our results from one sampling coincide with the exact 
expectation values. In section [4] we investigate the model on the 36-site lattice, 
for which a number of low-lying energy eigenvalues are known from the exact 
diagonalizations[2j. Starting with approximate states for the ground state of 
the system, which we obtain through procedures described in Appendix B, we 
argue the accuracy of the constrained SSS method by comparing our results with 
the ground state energy reported in ref.[25. The resultant expectation values 
show that the constraints are effective to improve sampling errors. Section [5] is 
to report our numerical results on the 48-site lattice where we know neither the 
exact eigenstate nor the exact eigenvalue. Assuming some symmetries which 
exist in the model on the 36-site lattice, we evaluate expectation values of the 
TO-th power of the operator Q obtained from several approximate states whose 
Sz, the z component of the total spin S of the system, is less than or equal to 4. 
We then estimate lowest energies for each sectors with S'z = k (k = 0, 1, 2, 3, 4). 
We see our data are well described by the arguments based on the spontaneous 
symmetry breaking. The final section is for summary and discussions. At the 
end of the paper we add three appendices in order to give detailed description for 
some parts of our numerical study. Appendix A is to explain how we determine 
values oi I in ^. In Appendix B we show our procedure to obtain approximate 
states on the 36-site and 48-site lattices. Finally Appendix C provides an em- 
pirical formula we use in the evaluation of the systematic error which is caused 
by employing the power method with a finite value of the power. 



2 Method 

In this section we describe our method. Subsection l2.1l is to give a brief descrip- 
tion of the conventional SSS method, the SSS method which is not accompanied 
with any constraints JJ^. Then subsection 12.21 foUows to show our basic idea for 
constraints. FinaUy in subsection 12.31 we describe the constrained SSS method 
in detail. 

2.1 The SSS method 

The stochastic state selection is realized by a number of stochastic variables. Let 
us expand a normalized state | ip) in an A^-dimensional vector space by a basis 
{I j)}j I V") = 12i=i I j)'^i- Then we generate a stochastic variable rjj following 
to the on-off probability function, 

P^{^) = ls{rj-a,) + il--)5{rj), 1 = min f 1, MV (3) 

aj aj aj \ e J 

A positive parameter e which is common to all Pjif]) {j = 1, 2, • • • , N) controls 
the reduction ratcQ Note that rjj — Oj (> 1) or ijj — and statistical averages 
are {{rij)) = 1 and ((77^)) ~ aj because of ([3|). A random choice operator M is 
then defined by 

N 
^^Elj>:'0-| ■ (4) 

Using this M we obtain a state | ip) , 

N 

which has less non-zero elements than | tp) has. We call the difference between 

I ^) and I V), 

|X)5^IV^>-|^), (6) 

a random state. 

Since {{rjj)) = 1 an expectation value {^p \ O \ ip) with an operator O is exactly 
equal to the statistical average (( {ip \ OM \ ip) )). 

Note that in the conventional SSS method all rjj 's are independently gener- 
ated stochastic variables. 

2.2 Basic ideas for constraints 

In this subsection we show our basic ideas in a very simple toy model. This 
model has only two basis states in the vector space, | 1) and | 2), which are 



^If the state | tjj) is not normalized, e in ||3) should be replaced by eyj {^ | tp) 



orthonormalized as (1 | 1) = (2 | 2) = 1 and (1 | 2) = 0. If we apply the 
conventional SSS method to a state | i/i) = | l)ci+ | 2)c2 (ci ^0, C2 7^ 0) in 
this vector space we obtain a state | ijj), 

I ^) =1 l)cir,i+ I 2)C2772, (7) 

where both 771 and 772 are stochastic variables generated by ([3]). We can repro- 
duce the state | V') by the averaging process since 

(( I ^) )) =\ l)cAm))+ I 2)c2((r;2» =| l)c,+ \ 2)c^ =| V> (8) 

Now we notify that independency between r]i and 772 are not needed in ([8]) 
because (( \ -ip) )) =| tp) is fulfilled as far as ((771)) = (irj2fl = 1- By making use 
of this possible dependency we can impose a constraint. For example, let 772 be 
not an independent stochastic variable hut the following function of 771 , 

^2 = l+(^) (l-^i). (9) 

It is clear that ([8]) holds with ([9]) because ((772)) = 1 follows from ((7/1)) = 1. 
With JD) we also see that 

(7/. I tA) - c?77i + cli^2 = cl + cl = {iP\^) (10) 

holds in each sampling. In other words, we have a constraint that a normaliza- 
tion {ij} I -0) is a constant c\ + c^. This also means that 

(^ I X)5 - (11) 

for any sampling, where | x)g is the random state defined by ([6|). 
It is possible to impose a more general constraint 

($ I X>5 = (12) 

instead of (HI]), with a given state | $) =| l)bi+ \ 2)62 (61 7^ 0, 62 7^ 0). To do 
so we should let 

Then we obtain 

($ I 7/i) = 61C1771 -I- 52C2r;2 = 61C1 + 62C2 = ($ I V-), (14) 

which is equivalent to p^ . 

2.3 The constrained SSS method 

In this subsection we present a way to impose constraints in the SSS method. It 
is straightforward to generalize discussions in subsection 12.21 for a larger vector 



space constructed by N basis states. Let the constraint be, with a given state 

($ I x)9 = (15) 

for a state | ip) = J2j=i I J)'^j ^^'^ I x)9 = S7=i I J)'^ji''lj ^ !)■ I^^ this case we 
should solve the equation 



TV 

E 



bjcjiv, -1) = (16) 



to give one of rjj , say ryj , as a function of A^ — 1 independent stochastic variables 
with j 7^ J. 

From the fact that {{rjj)) = 1 for all j except for J, it is clear that pT]) guarantees 
((7?j)) = 1. As for ((772)), we find 

((.^»-1 = E(^) '(%-!)• (18) 

In principle, J can be any of 1, 2, • • •, N. From practical point of view, however, 
J should be chosen carefully so that sampling errors are diminished. When we 
choose the state | ip) as the state | $), which means bj = Cj for all j, the right- 
hand side of p8)) becomes "^j-ij ( ^ ) (oj — 1)- In this case it is clear that we 
can lessen this quantity by picking up J which realizes \cj\ > \cj\ for any j. 
Numerical examinations for this choice will be given in following sections. 

Another generalization to impose more than one constraint is also easy. Let 
K denote the number of constraints. As far as if < iV we can impose constraints 
for given states | <E>fe) {k = 1,2,- ■ ■ ,K) 

($fclx>5 = 0, (fc = l,2,...,X), (19) 

by requesting K variables among rjj^s depend on N — K other rjj^s which are 
independent stochastic variables generated by ([3|) . More concrete description of 
the way to impose constraints is as follows. If we know the coefficients of the 
expansions of these states, the above constraints (fTO)) read 

J2bfc,{vj-1) = {k = l,2,---,K), (20) 

3 

with 

N N 

i^> = Ei^>^-' ii'fc> = Ei^>f- (21) 

From ([20]) we obtain 77,7^., (fc — 1, .., K) which are dependent on other stochastic 
variables. Since these dependencies are linear, the conditions ((77,7^,)) = 1 are 



always satisfied. A choice of {Jk} is arbitrary, but the sample fluctuations will 
depend on the choice. 

Hereafter we denote a random choice operator used in the constrained SSS 
method by Mc in order to avoid confusions with the random choice operator in 
the conventional SSS method. 

3 Study with an exact eigenstate 

In this section we study a system for which we know an exact energy eigenstate 
as well as its eigenvalue. First we make analytical discussions with the exact 
eigenstate | iPe) whose exact eigenvalue is E. Then we numerically study the 
spin 1/2 quantum Heisenberg antiferromagnet on a 16-site triangular lattice 
as a concrete example. For this small lattice we can easily obtain the exact 
ground state and its eigenvalue by the exact diagonalization. This enables us 
to compare our results from the constrained SSS method with the exact ones. 
Let us examine the expectation value of the vn-th power of Q = II — H, 
(V' I Q™ I ip), with I "0) =1 V-'e)- From the exact eigenvalue we obtain 

(V-E I g" I ^e) = Q"\ Q = l-E. (22) 

In order to calculate these expectation values by the constrained SSS method, 
we insert random choice operators Mc (n = 1, 2, • • • , m) and calculate 

(^EigMi^'-.-QA/WlVE). (23) 

Here we denote random choice operators by Mc instead of Mc since we want to 
emphasize that different operators, each of which includes stochastic variables 
independent of those in other operators, are used. Let us define states | i/)'^"^) 
as follows, 

|0(")) = gM<("-i)...gMW|V'E) (n>2), |0W)=|Ve). (24) 

Note that I (/)(")) = gAfi""^^ | 0("-i)) with n > 2 by definition. For each Afi"^ 
we impose the dependency (J17l) so that the constraint 

(V-E I X^"^)5^"^ = (25) 

holds for 

I x("))g(") ^ Mi") I 0("V I 0^"))- (26) 

Here J*^"^ of the dependent stochastic variable is determined by the condition 

|cV(„)| = maxi<j<jv |Cj" I with the expansion | 0(")) — ^ | j)Cj . Notifying 
that (V'E I M""^ I 0^")) = (V'E I <^^")) follows from ^, we obtain 

(V'E I gA^i") I </>«> = giV-E I M(") I 0(")) = g(VE I <^("') 

= g(V.E I gA/i"-^) I 0("-^)). (27) 



Using (P7)) repeatedly, we find 

(^E I QMi"^ • • • QMi'^ I ^e) = (V'E I QM^"^^ I 0^")) 

- Q(V^E I QMi"-^) I 0^"-^)) 

= Q^-'iV^E I Qm(i) I 0(1)) 

- Q™(Ve|mW|0W)=Q™(V'eU(^^> 

= Q'^iV'E I ^e) = Q". (28) 

This result implies that the exact calculation is possible in the sampling. Note 
that only one sampling is enough in the constrained SSS method. 

Now we present numerical results for the ground state of the system on the 
16-site triangular lattice. The ground state energy is known to be £■ = i?g = 
—8.5555. Since we employ I — 2 here as is stated in Appendix A, the exact 
value oiQ \sl — Eg = 10.5555. In Fig. [T] we plot ratios of the expectation values 
(ED, 

„(,„) ^ (^E I 0<"+i)) ^ (^E I QMJ"'^ I 0'")) .2q. 

- (V'eU^™)) (^eU(™)) ' ^ ' 

obtained in one sampling with the constrained SSS method. The parameter 
e in ([3]) is 0.05. As is expected from above discussion, we observe i?*^™) = Q 
for any value of m. We also present results obtained by the conventional SSS 
method, averages of (?Ae | QAf (™) • • • QM^i) | V:e)/(V'e | QAf (™-i) • • • QAf (i) | 
V'e) from 100 samples with e — 0.05, where M^^^'s denote different random 
choice operators in the conventional SSS method. Statistical errors for these 
averages are also plotted in the figure, where for data Xi (i = 1, 2, • • • , rigampic) 
the statistical error A[A'] is defined by 

A[X] = ^ '^[^J , (30) 



Y' ^sample ^ 



with the standard deviation 



a[X]^ 



faampjc fcampic 



\ ^sample _, ^ ^ 



sample 



(31) 



Comparing these data we clearly see that fluctuations existing in the conven- 
tional SSS method disappear with constraints stated by (PS)) . 

4 Study with an approximate eigenstate 

In this section we argue the case in which we use an approximate eigenstate 
I V'a) instead of | V'e)- First we present an analytical argument which endorses 
that the constrained SSS method is effective for eigenvalue evaluations starting 
from an approximate state. Then using the 36-site lattice we demonstrate that, 



starting with an approximate ground state, we can obtain in the constrained 
SSS method the correct ground state energy with much less fluctuations than 
those from the conventional SSS method. 

Let us start with an approximate state | V'a) which has some overlap with 
the corresponding exact eigenstate | i/^e), 

\M=w\^e) + s\0- (32) 

Here | ■0a) is normalized with w'^ + s^ = 1 and we expect \w\ ^ \s\. Instead 
of (P5|) we calculate the expectation value {tpA \ QMc ■ ■ ■ QMc \ V'a) with 
constraints 

(0A I M^'^^ I 0i"^) = (0A I 0i"^> (n = 1, 2, . . . , m) (33) 

where 

|0i"))^QAfi"-i)...QMW|0A) (n>2), | ^i^^) ^| 0a). (34) 

Let us examine ('0a | QMc \ 0a ) when the constraint (|33p holds. Using (|32p 
we find that 

(0A|QAfW|0W) ^ (u,(^E|+s(C|)gMW|</'i'^) 

= (q«;(^e I +s(C I Q) M(i) I 0i")) 

= (g[(0A|-.(CI]+s(CIQ)MW|0«) 

= g(VA I M^^^ 1 0i'^) + s(c I (q - gi) M^ 1 0i'^) 

= g(V'A 1 0i'^) + s(c I (g - Qi) (i 0i' V I xi'^)5i'^) 

= g(V'A I v-a) + s(c I (g - Qi) (i ^a)+ I xi'^)5i'^) 

= g + s(ci (g - g/) i ^a) 



S{ 



.W\.W 



(CI g-g/ lxrH'^ (35) 



where 



ix!r^)gr^M(")|0i"Vi0i"^ (36) 

Note that in the right-hand side of ([5S)l only the last term contains the fluctu- 
ation by a sampling and that this term should be small when \s\ <^ 1. 

Next we turn to the numerical study for the ground state on the 36-site 
lattice. The exact ground state energy of this system is known to be —0.186791 
per bond_24_, namely —20.1734 in total. Under the symmetries the ground 
state of this system has, the number of basis states in the whole 3^=0 sector 
amounts to ~ 2.2 x 10^. Following the procedures described in Appendix B 
we create two | 0'a)'s for the ground state of the system. Let us denote them 
by I 0A/i) (m — 1^2). The number of basis states with non-zero coefficients 



in the expansion of | V'Ai) is 333001 and the expectation value of H is (V'ai | 
H I Vai) = —19.577, while | ?/'A2) includes 887875 basis states with non-zero 
coefficients and ('0A2 | H \ ipA2) = —19.817. In the same manner as was stated 
in section [31 we request one constraint for each random choice operator. We 
choose J^"^ for the only non-independent variable rij(n) using the same criteria 
as that in the Ns ~ 16 case. Based on conditions we notify in Appendix A, 
we choose I = 4 (namely Q — il — H) here. Our results for the system are 
presented in Figs. [l[3ll|and[5l Fig. H plots 

^(™) ^ (V'A I <t^'r"^) ^ (^A I QMi"' I 4>t^) 

(^Al^r^) (^A|0i™^) 

calculated from 100 samples with the constraint stated by f33)) . The value of the 
parameter e is 0.01. In the figure we also plot values of (V'a | QM^™"* ■ ■ ■ QM'^^^ \ 
V'a)/(V'a I Q-^^™"^' • ■ •QAf(^) I '0a)i which are obtained by the conventional 
SSS method from averages of 100 samples with e = 0.01 using the approximate 
state I V'Ai)- We see that the values mostly agree in both method. In Fig. [3] we 
compare the standard deviations a of the ratios shown in Fig. [5] We see the 
exponential growth of the standard deviations except for the m = 1 datum from 
the constrained SSS method. It is quite impressive that fluctuations from the 
constrained SSS method are much less than those from the conventional SSS 
method. 

Let us present some more data on the 36-site lattice which will be helpful to 
understand the role of the parameter e and the quality of the approximate state 
I ^a) • Fig- HI is to show how much basis states are included in the expansion of 
Mc I 0a)- Here we present only results from | V'Ai) since we observe A'a's 
are mostly insensible to the choice of the | 0a ) =1 ^a) stated above. We see 
that as m increases iVa becomes almost constant for each value of e and there 

iVa ex e-i (38) 

holds. Note that this means that by the choice of e we can change the CPU 
time and the memory which we should supply in numerical studies. Fig.[5]shows 
several values of Rj^ ' from 10 samples obtained with different values of e as 
well as with two different approximate states | V'Ai) and | ip a2)- It is quite 
reasonable that the data indicate we can obtain a better lower bound for Q (a 
better upper bound for E) with a better approximate state | V'A2) when the 
value of rn-max is the same. In order to prepare a better approximate state, 
however, we of course have to deal with a larger portion of the Hilbcrt space. 
We also see in Fig.[5]that statistical errors, and therefore the standard deviation 
which is 3A[_R^ ] in this figure with Y^nsampio — 1 = 3, are irrelevant to the 
choice of the approximate state, while they decrease rapidly for less value of e. 
We observe that 

a[R^l'^]^e'< (7-3). (39) 
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5 Application to the 48-site system 

In this section we study the ground state energy and some excited energies 
of the spin 1/2 quantum Heisenberg antiferromagnet on the 48-site triangular 
lattice using the constrained SSS method. Similarly in sections [3] and 31 we 
request one constraint for each random choice operator and decide J^"^ for the 
only non-independent variable ?7j(n) with the condition |c,"„)| > \cl \ for all j. 
In each sector with S^ — n [n — Q^ 1,2,3,4) we study the state whose energy 
eigenvalue is the lowest in that sector. Following the procedures in Appendix 
B we calculate | '(/'a)'s summarized in Table 1. As is stated in Appendix B, 
symmetries we assume for the approximate states are the same as those for the 
Ns — 36 case. Figure |6] shows our results for the ratio R™ defined by p7)) 
up to TOinax — 15 with valucs of I given in Appendix A. All data are calculated 
from 100 samples with e = 5 x lO'^ (e = 7 x lO'^ ) for 5^ = (5*^ > 0). 

It is known that the argument by the spontaneous symmetry breaking in 
semi-classical Neel ordered antiferromagnets suggests the following energy spec- 
trum on finite-sized lattices [3l l28l [34l [35] . 

E{S)-E{0) = -^-^±^, (40) 

2XA Ns 

where E{S) and xa denote the lowest energy of the system with the total spin 
S and the susceptibility, as long as S* ^ \/Ns. Keeping this in mind, we plot 
values for Ef^/{3Ns) versus S{S + 1) in Fig. [7] with an assumption S ~ k, where 
Ek. denotes the upper bound of E with Sz ~ n given by /j^"'"'*^-' . Namely, £« 
is I — R)^ in the Sz = k sector. Values of iJ^ are also presented in Table 1. 



We see that those energies on the 48-site lattice are well described by (|40|) . 
From the least square fit of the data with S — k = 0, 1, 2, 3 and 4 we obtain 
1/2XA = 5.6. In Figl7]we also plot the data on the 36-site lattice obtained by 
the exact diagonalization[24] . which gives 1/2xa = 5.1. The fact that this value 
is compatible with the discussion for the susceptibility in ref. [28] supports the 
finite size arguments ()40p based on the symmetry breaking. 

Finally let us comment on the spin gap Aspin(A^s)- Results on Eo and Ei 
give us Aspin(48) — 0.284 ± 0.027. Through the finite-size extrapolation using 
data Aspin(12) Aspin(36) and Aspin(48) we obtain Aspin(oo) ~ 0.10, which is 
smaller than the one evaluated from the data with Ns < 36 [3|. 
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s. 


Ma 


Ma/N^ 


(V-A 1 i? 1 V-a) 


e 


i?. 





4.9 X 10« 


8.8 X 10"^ 


-26.411 


5 X 10-^ 


-26.611 ±0.015 


1 


1.5 X 10» 


2.8 X 10-^ 


-25.981 


7 X 10"^ 


-26. 327 ±0.023 


2 


1.2 X 10« 


2.4 X 10-^ 


-25.482 


7 X 10-^ 


-25.815 ±0.020 


3 


7.3 X 10^ 


1.9 X 10-^ 


-24.753 


7 X 10-^ 


-25. 064± 0.017 


4 


9.0 X 10^ 


3.1 X 10-^ 


-24.094 


7 X 10-'^ 


-24.302 ±0.016 



Table 1: Approximate states employed for the 48-site system and results for the 
ratio i?^""""-' with mmax = 15, which we obtain from 100 samples using the 
constraint ([33]). A^„ denotes the number of basis states of the whole Sz — k 
sector with the assumed symmetries, while A/a is the number of basis states 
whose coefficients are non-zero in the expansion of | V'a)- 



6 Summary and discussions 

In previous works [T51 [TOl [Tfl [TBI UM \7U\ we have developed the stochastic state 
selection (SSS) method, which is a kind of Monte Carlo method suitable to 
calculate eigenvalues in large quantum systems. In this paper we proposed 
the constrained stochastic state selection method, where some constraints are 
imposed in each sampling. 

It is a characteristic feature of the SSS method to sample many states si- 
multaneously. Namely, in the SSS method we introduce N stochastic variables 
to form a random choice operator, where N denotes the number of the basis 
states of the whole vector space. In the conventional SSS method all of these 
N stochastic variables are independent of each other. In the constrained SSS 
method, on the other hand, K variables in the random choice operator are de- 
termined by values oi N — K other ones which are independently generated 
stochastic variables. Using these dependencies we can force K relations to rep- 
resent constraints should hold in each sampling, provided that these relations 
are linear of all stochastic variables. 

We pursue arguments to a conclusion that the constrained SSS method is 
effective to decrease statistical errors in calculations of energy eigenvalues using 
approximate eigenstates. Numerical demonstrations then follow, which apply 
the constrained SSS method to the spin-1/2 quantum Heisenberg antiferromag- 
net on a 36-site triangular lattice. Here we employ the power method in com- 
bination with the constrained SSS samplings. We impose one constraint which 
requests that each sample for a given state should not change the inner product 
between the initial approximate state and the given state. With initial states 
which approximate the ground state of the system, we calculate expectation 
values to be used in the power method. We observe much less fluctuations in 
the constrained SSS method compared to those in the conventional SSS method 
and the ground state energy estimated in the constrained SSS method are in 
good agreement with the known exact value. 

We then successfully calculate the lowest energies in Sz = k {k = 0, 1, 2, 3, 4) 
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sectors on the 48-site triangular lattice. Our results on low-lying energies with 
different values of Sz add an evidence to the ordered ground state and the 
finite-size arguments based on the symmetry breaking. Especially, our result on 
the spin susceptibility obtained from the ground state energy and the low-lying 
excited state energies is consistent with finite-size effects reported in ref. pS] . 

Several comments are in order. 

In study of the Sz — sector in ref. [19] we use the Lanczos method together 
with the conventional SSS method. In the present work we employ a simple 
power method instead of the Lanczos method so that the calculation algorithm 
is simple enough to enable us to calculate eigenvalues with Sz =/= within our 
computer resources. 

The efficiency of the method in numerical studies is mainly controlled by the 
value of the parameter e. In general we can expect better results with smaller 
values of e, but those calculation would take more CPU time and memory. 
Therefore in actual calculations one should choose the value of e so that the 
computer resources stay within limits of his computerf|j. Roughly speaking, the 
memory size Mcpu is proportional to TVb , the number of basis states which have 
non-zero coefficients in the expansion of | 0'™+^') = QMc \ '/)''"■'). Since we 
observe that A^b is proportional to A^a shown in Fig. [4] and that its dependency 
on e is described as (155)) . it leads 

Mcpu oc e-\ (41) 

The total CPU time Tcpu, on the other hand, has a more complicated relation 
with e because it depends on the number of samples ^sample as well as Ng_. It 
also depends on the number of iterations rrimax when we use the method in 
combination with the power method. Therefore 

TcPU OC m,„ax • ?^samplo ' Na- (42) 

Values of m,„ax and risampic are determined from the results for the standard 
deviations a such as those shown in Fig. [3l Since, as we have observed an 
exponential growth of cr as a function of m (except for to = 1) there, ffuctuations 
of the data grow rapidly when iterations are repeated. We therefore have to 
give up our numerical study with some finite value of m before the data become 
statistically meaningless. If we can employ a smaller value of e, the value of 
'm-ma.x would become larger because of much improvement of a suggested by 
([39]) . The number of samples should be chosen so that the statistical errors are 
reasonably small for m < TO,max- The total CPU time we spent to calculate 
the data presented in this paper is about 2000 hours with a computer whose 
memory is 8 Giga Bytes and whose CPU is Xeon Dual Core. 

What can we say about the accuracy of our results? Let us here estimate 
the systematic error which exists owing to the power method with finite powers 
up to TOjnax- On the 36-site lattice the exact ground state energy is known to 



^ Note that the method is appHcable even with large values of e provided that the number 
of samples is large enough. 
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be Eg = —20.1734. The upper bound of the systematic error therefore can be 
estimated by {Eu\y — Eg)/\Eg\ with an £^ub, an upper bound of i? given by l — Qih 
where Qib denotes the lower bound of Q. We see that _R^°"'='-' gives our best 
lower bound of Q, because R^ increases as m grows and R^ — > Q {m — > cx)) 
should hold. Using I — Rj^ ' in Fig. [5] as E^h we conclude that the systematic 
error is less than 0.8% (0.5%) for | V'Ai) (| V'A2)). On the 48-site lattice where 
no exact energy eigenvalue is known, we try to find a lower bound Ei^ from 
our data using the upper bound E^^ of E. Note that this task is equivalent to 
find Qub, an upper bound of Q, because E\iy ~ I — Quh- In order to find an 
upper bound of Q we carry out an empirical fit. Details of this fit are described 
in Appendix C. Then using Eub and Eih we evaluate the systematic error by 
{Euh — E\h) /|£^ub| in each sector with Sz = n- Results from the data presented 
in Fig. [3 with Qib = i?^^""""') = rf^) are 0.6%, 0.9%, 0.8%, 3.2% and 1.3% for 
Sz = 0, 1, 2, 3 and 4, respectiveljH 

It would be the simplest way to impose only one constraint in each random 
choice operator {K = 1) as we did in sections [3] and H) Nevertheless one should 
remember that imposing more constraints in one random choice operator {K > 
1) is also, at least theoretically, possible in the constrained SSS method as was 
discussed in section [51 Although the A' > 1 calculations might be numerically 
more difhcult, further study for these cases is desired from practical point of 
view. 

Finally let us emphasize that the constrained SSS method, as well as the con- 
ventional SSS method, has no physical bias since this method does not depend 
on any physical assumption Q Therefore the method is applicable to numerical 
study of various systems. Results obtained in this work on triangular lattices 
encourage us to numerically study spin systems on other lattices by means of 
this method. 
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Appendix A 

Here we describe how we determine the value of I in the operator Q defined 
by dH). Let us denote all eigenvalues of iJ by ii^min(< 0), E^, E^,, ■ ■ ■, -Emax(> 0), 
where 

^min < Ea < Eh < ■ ■ ■ < i^max- 



^ Similar analysis on the 36-site lattice yields the result that, with mmax = 20, 
{Eub - Eli,) /\Eub\ is 0.9% (0.4%) with | Vai> (I V'A2)). 

■^Although in applications presented in this paper we assumed some symmetries to construct 
approximate states, these assumptions are not essential in the constrained SSS method itself. 
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It is easy to see E'max = 3iVs/4 for the TV^-sitc lattice. Since Q = II — H, 
eigenvalues of Q are then 

I — £^min >l — Eg,>l — E\,> ■ ■ ■ >l — i?max- 

Let us consider a state | VP) and expand it using an orthonormal basis { | 'I'x) } 
which is obtained from eigenfunctions of H , 

I *> = E I *-)/- 

X 

Then we obtain 



2/7 77 ^™ n I Y^ 



f^ii-E,r{i + 



fl fl-E. 



t^f^y^-^y 



f^ii~E,r{i+j2ll(i-''- ^^ 



x/y 



f^\ l-Ey 



where by suffix y we denote the term whose |Z — i?x| is the largest among x with 
non-zero /x. It is clear that the term /^ {I — -Ey)™ dominates as m increases. 

When I *) =1 ■0a) defined by ^^ with f^i^ \ *min) = w \ V'e), the value of 
I therefore should satisfy the condition 

1^ — -EminI > 1^ ~ -EmaxI 

in order for us to pick up the term with f^^^^^ {I — -Emin)™ in (V'A I Q™ I V'a)- 
Limiting ourselves to the range I < -Emax we thus see / + |E,„in| > En^ax ~ I, 
namely 

I > X (Emax — |Emiii|) 

should hold. In most cases we have to use an upper bound of E,„i„, Euppcr, 
instead of Emin. Then we determine the value of I in the range 

'' -^ 7\ l-^niax |-^uppcr|j i 

which includes the range I > ^ (E^ax ^ |Emin|) for any Euppor < because 

Emax — I Euppcr I > E^ax " |Emin| holds. 

In choosing a value of I which satisfy the above condition, one should also 
note that contributions from excited states increase as / increases because the 
dumping factors ( 1 — ^1"^°'" 1 decrease. 

For the 16-site lattice, we use the exact value Emin — —8.5555 to decide I = 2. 
For 36-site and 48-site lattices we use Euppcr calculated from (0a | H \ 0a)- 
Values of / we chose are summarized in the Table 2. 
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Ns 


Sz 


-C/uppcr 


-^inax |-^uppcr| 


/ 


36 





-19.57 


7.43 


4.0 


-19.81 


7.19 


4.0 


48 





-26.41 


9.59 


5.0 


1 


-25.98 


10.02 


5.2 


2 


-25.48 


10.52 


5.5 


3 


-24.75 


11.25 


5.8 


4 


-24.09 


11.91 


6.2 



Table 2: Values of / we use in the Sz = k sectors {k — for the 36-site lattice 
and < K < 4 for the 48-site lattice), where Sz is the z component of the total 
spin jS*. -Euppcr is an upper bound for the lowest energy eigenvalue of the system, 
which we obtain from the expectation value of the Hamiltonian calculated with 
an approximate state, while i^max denotes the maximum energy eigenvalue of 
the system. 



Appendix B 

Here we explain how we obtain approximate states | '0a)'s on the 36-site 
and 48-site lattices. The method is essentially the same as the one given in 
ref. [19, for the 36-site triangular lattice. The only difference is that we include 
here as many degenerate Ising-like configurations as possible in an initial trial 
state. This improvement comes from Wannier's rigorous proof[36j which claims 
that a classical antiferromagnetic Ising system, namely the spin system at zero 
temperature T = 0, on an A^s-site triangular lattice is heavily degenerated for 
its minimum energy —Ns/i. 

For our numerical work in this paper, we employ one basis on the 36-site 
lattice and five bases on the 48-site lattice corresponding to values of Sz ■ First 
we comment on these bases together with brief descriptions for the transforma- 
tion symmetries we impose on them. Then we show procedures to create an 
approximate state for numerical studies. 

Each state | j) in our basis states {| j)} is represented by linear combinations 
of states I si,S2, ■ • • , satJ with s„ = -1-1/2 or —1/2 (n = 1, 2, • • • , Ng), which 
indicates that the z-component of the spin on the n-th site is -t-1/2 or —1/2, 
respectively. They belong to one of 2Ns + 1 sectors according to the value of 
Sz = X]n=i ^n, which is the z-component of the total spin S of the system. 

The basis states to form the approximate state for the ground state on the 
36-site lattice should have symmetries found in ref. [24]. Following them we 
construct our basis states with 5*2 = so that they have translational symme- 
tries for the zero momentum (0, 0) and even under the 27r/3 rotation and the 
reflection. Each space group which consists of the translation group and the 
point group has then 432(= 12 x 36) elements. 

For the model on the 48-site lattice we assume each basis state with a given 
value of Sz has the same symmetries as those of the exact state for the Ng = 36 
lattice [24] whose energy eigenvalue is the lowest in the corresponding sector. 
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Namely, all basis states in both bases have translational symmetries for the 
zero momentum (0, 0) and even under the ir rotation, the 27r/3 rotation and the 
reflection. The basis states in the basis used for S'^ = 0, 2 and 4 calculations also 
have the even tt rotational symmetry, while those in another basis prepared for 
5*2 = 1 and 3 have the odd one. Both space groups have 12Ns = 576 elements. 

Now we start to calculate | ■(/'a) in each sector with a fixed value of Sz, using 
an appropriate basis among those stated above. 

The first stage is to find as many degenerate states for the eigenvalue —Ng/A 
as possible. Using the conventional Monte Carlo method at low temperature(r = 
0.5), where the classical energy is used as the Boltzmann weight, we pick up 
iVit -- 10* states to fulfil the condition {j \ H \ j) = -N^/i and S^ = n 
(0 < K < 4), where the Hamiltonian H is given by ([T]). With this Na basis 
states we calculate the eigenstate | \E't) for the lowest energy eigenvalue within 
this partial Hilbert space by means of the conventional exact diagonalization. 

The next stage to calculate an approximate state | ipA) is to repeat following 
procedures, starting from the initial | \l/t) obtained above, until the expectation 
value for H does not change beyond our criteria. 

i. Extend the partial Hilbert space given by | ^t) through operations of a 
few H's to I ^t), until the available computer memory is exhausted. The 
maximum number of basis states we can permit is about 10^. 

ii. Within the Hilbert space determined in (i), pursue the state | ^t') with 
which (^t' I H I 'I't/) is as low as possible. 

iii. If the change of the obtained value (^t' | H \ ^t') is in the range of five 
decimal digits, employ | Vl/f) as the approximate state | V'a)- Otherwise, 
form a state | *t") by keeping a few percent of the basis states whose 
coefficients in the expansion of | ^t') sue relatively large and replace | ^t) 
by I ^t") to proceed with (i), (ii) and (iii) once more. 



Appendix C 

Here we explain how we evaluate the lower bound of the eigenvalue of Q = 
II — H, which we denote by Q\h, from the data Rj^ {ni = 1, 2, • • • , mj^ax) and 
a given value of the upper bound of the eigenvalue of Q denoted by Quh- We 
employ an empirical formula we presented in a previous paper [15j, which is 

F(m,Qw,go,a) = Qw fgo+ J\ , ] ^ qi = {a + 1){1 - qo). (43) 
V m + a + 1 J 



The quantity we measure is •\/D(77i,„ax, Qw, qo, ot), 

(^A 1 0!r^) 



D{m,^^^, Q^, qo, a) = ^ 



771 — 1 
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F{m,Qy,, go, a) 



(44) 



where values of {ipA \ 4>a ) ("^ = Ij 2, • • • , rn-max) are calculated from the data 

771—1 

{^A\cbt^)^Y[R^:\ <^1. (45) 

Changing Qw around Qib which is given by the data i?j^""""' , we look for val- 
ues of parameters qo and a, say q and a, so that ^(TOmaxiQw, 9„j a) gives 
a local minimum for the given value of Qw Fig. [5] shows a typical result for 
< /D(m,„ax, Qw, 9r|)fi) as a function of Qw By requesting that the value of \/D 
at Qw = Quh should be equal to the value at Qw = Qib, we obtain a upper 
bound Quh shown in the figure. 
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Figure 1: Ratios defined by ([29]) calculated on a 16-site triangular lattice 
with e = 0.05. Here we use the exact ground state as | Ve) for which 
E = Eg = —8.5555. The crosses plot data from one sample obtained 
by the constrained SSS method. We also plot, by open diamonds, ratios 
(V'E I QMm- ■ -QMi I Ve)/(V'e I QMrn-1- ■ -QMi \ Tp^) which are averages 
from 100 samples generated by the conventional SSS method. Errors in the 
figure are statistical errors only. 
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Figure 2: Ratios ([37]) on a 36-site triangular lattice which are generated by the 
constrained SSS method (filled marks). We employ two states to approximate 
the ground state energy. By filled circles we plot results which we obtain using an 
approximate state constructed by 333001 basis states (| ■0Ai))- With the same 
approximate state we also calculate ratios (V^a | QM^"^> ■ ■ ■ QM^^> \ iIja}/{iPa \ 
QM(™-i) . . .QM(i) I ^a) from 100 samples by the conventional SSS method 
with e = 0.01, which we show by open circles. Another approximate state 
constructed by 887875 basis states (| ipA2)) is also used in the constrained SSS 
method, whose results from 100 samples with e — 0.01 are plotted by filled 
squares. The dotted line in the figure shows the exact value of Q obtained from 
the exact ground state energy [24]. Errors in the figure are statistical errors only. 
Statistical errors for both filled marks are within marks. 
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Figure 3: Standard deviations for ratios shown in Fig. [21 which we obtain using 
the approximate state constructed by 333001 basis states (| V'Ai)), from 100 
samples with e = 0.01. Filled circles present a[R™ ] obtained by the constrained 
SSS method. Open circles are standard deviations for ratios obtained by the 
conventional SSS method. 
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Figure 4: TVq, the number of basis states with non-zero coefficients in the 
expansion of Mc | 4>k )■ ^^1 data are calculated from one sample with 
I ipfi^ —\ tpAi)- Pluses, diamonds, triangles and asterisks present the results 
with e = 5 X 10^^, 1 x 10~^, 5 x 10^^ and 1 x 10^^, respectively. 
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Figure 5: Ratios /j^'"""'-' with vTimax = 20 versus e. Circles (squares) present 
the resuhs obtained from 10 samples with | V'Ai) (| V'A2))- Errors shown in 
the figure are statistical errors only. The standard deviations calculated from 
these 10 samples are consistent with those shown in Fig. [5J where the number 
of samples is 100. The dotted line indicates the exact value of Q obtained from 
the exact ground state energy [24]. 
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Figure 6: Ratios ([37]) in the Sz — k {0 < k < 4) sectors on a 48-site triangular 
lattice calculated from 100 samples with e = 5 x 10^"^ (e = 7 x lO^'^) for 
Sz = {Sz > 1). We use approximate states which are composed of 489413140, 
150733425, 115759910, 73294432 and 90008649 basis states in the Sz = 0, 1, 2, 3 
and 4 sectors, respectively. Circles present the Sz = results, while triangles- 
up, -down, -left and -right show results for Sz — 1,2,3 and 4. All statistical 
errors are within the marks. 
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Figure 7: Energy per bond obtained by Ek/{SNs) {k = 0, 1,2,3,4) on the 48- 
site lattice. The data are shown by asterisks, whose errors are within the marks. 
We also plot the 36-site data from ref. [53] by crosses. The dotted (dashed) line 
is obtained by the least square fit assuming ([40]) on the 48-site (36-site) lattice. 



24 




Qw 



Figure 8: Typical plots of jD{mn^s^^,Qy,,q^,a), where D(mmax, Qw, go, a) is 

defined by ([33]), as a function of Qw Using this figure we can find Qub, which 
is indicated by a filled diamond, from a given value of Qib indicated by an open 
diamond. The dotted horizontal line is to guide eyes. 
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